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We study states with spontaneous spin current, emerging in frustrated antiferromagnetic spin-S 1 chains subject 
to a strong external magnetic field. As a numerical tool, we use a non-Abelian symmetry realization of the 
density matrix renormalization group. The field dependence of the order parameter and the critical exponents 
are presented for zigzag chains with S = 1/2, 1, 3/2, and 2. 
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1. Introduction 

In usual magnets, spin degrees of freedom order ferro- or antiferromagnetically, while those 
two basic orders may be viewed as special cases of a general spin density wave order parameter. 
It had been known for a long time [ [TJ [2] that other "exotic" spin orderings with a vanishing 
expectation value of the spin density are theoretically possible, but only recently such orders have 
been observed in numerical studies of realistic models [ [31 SI El [HI [Z] • One prominent example of 
unconventional order parameter is the spin current (the so-called vector chirality [ H] , or p-type 
nematic [[2] order parameter). In a spin chain, the vector chirality operator k(ji) = K,i, n +i can be 
defined as a vector product of two adjacent spins: 

4,n+l = {S n X S'n+l), (1) 

where n labels the lattice sites along the chain. Obviously, the vector chiral order is nonzero in 
a state which has a classical helical magnetic order, and the projection of k onto the helix axis 
distinguishes between the left and right spirals. For low-dimensional magnets, such unconventional 
order parameters may become dominating, since quantum fluctuations have a tendency to destroy 
the usual magnetic order. A state with vector chiral order may be viewed as the result of phase 
fluctuations destroying the helical order, but still keeping the preferred sense of rotation contained 
in a classical spiral [ [8] . The interpretation of k as the spin current can be understood from the 
equation of motion, which for a Hamiltonian of the form H = J2nm Jm,S n 1 S n+m reads Ti{dS n /dt) = 
,n+m}- One can show that in the ground state the net spin current should be 
zero, which for the above Hamiltonian translates into the condition yV mJ rn K n ^ n+rn = 0. 

Spin current states have been relatively well studied for anisotropic frustrated chains [ [3] 
[Till 21 [TT] . recently, it has been shown [ [T^l G3 [T31 H] that a strong external magnetic field can 
be used as a control parameter to drive the spin current in isotropic frustrated spin chains. The 
interest to this topic is further boosted by the fact that there are currently several quasi-one- 
dimcnsional magnetic materials that are intensely studied as possible candidates for manifestation 
of unconventional orders, particularly LiCuV04 [03], Li 2 ZrCu04 [DH], and Cu2Ci4-H 8 C4S0 2 (Sul- 
CU2CI4) [US]. We have previously studied [[3] the spin current correlations in frustrated chains 
with the spin 5=1/2 and 1, and have shown that the behavior of S = 1/2 and S = 1 chains 
was very different. It was not clear whether this can be attributed to some general difference in 
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the behavior of chains with integer and half-integer S. The goal of the present paper is to answer 
this question. For that purpose, we numerically study the spin current in isotropic spin-5 zigzag 
chains with the spin up to S = 2. It is shown that chains with S > 1 behave similarly and exhibit 
a state with a finite spontaneous spin current in the entire region of finite magnetization. In Sect. 
[2]we give a brief overview of the theory of spin current correlations, Sect. [3]presents the numerical 
results, and Sect. [4] contains a short summary. 

2. Frustrated chain in magnetic field: spin current correlations 

We consider a frustrated antiferromagnetic spin chain described by the following Hamiltonian: 
H = Ji S n ■ S n +i + J2 S n ■ 5„+2 — H S z (2) 

n n n 

where S n are spin-5* operators at the n-th site, J\ > and J2 > are the nearest and next-nearest 
neighbor exchange constants, respectively, and H is the external magnetic field, assumed to be 
applied along the z axis. At zero field, if the frustration parameter a = J2/J1 is large enough, 
the ground state is generally expected to have a finite spectral gap A, both for integer and half- 
integer S. When the applied field exceeds the critical value H c — A, the system acquires finite 
magnetization. Further increase of the field beyond the saturation field value H s brings the system 
into a fully polarized state. 

Let us first assume that a ^> 1, so that one deals with the limit of two weakly coupled spin-S" 
subchains. The behavior of spin chains in magnetic field have been extensively studied for 5* — 1/2 
[ [HI [H QI] and S = 1 [[201 GEO HI]- We will assume that for H c < H < H s the low-energy physics 
of a single chain is well described in terms of the effective Tomonaga-Luttinger liquid (TLL) theory 
with the following Hamiltonian: 

Htl[8,<p]=1 J dx[^{d xV f +K{d x ef). (3) 

Here K is the so-called TLL parameter, v cx J2 is the Fermi velocity, <fi is the compact bosonic 
field (ip = ip + v^r), and 9 is its dual satisfying the commutation relations \^p{x), 9(y)] = iQ(y — x), 
where Q(x) is the Heaviside function. 

We assume further, by the analogy with the S — 1/2 and S — 1 case [ |T7l [20j [22] , that in the 
continuum limit the most relevant low-energy part of spin operators can be represented as 

2 

S z a {x a ) = MS + ^dxipaixa) +A a sm[k F x a + V^K(pa(x a )] (4) 
V 77 

S+(x a ) = e™ x / 2 e^ e ^{A 1 +A 2 sm[k F x a + V4^ a (.T a )] }, 

Here a = 1, 2 labels the two subchains, the space coordinate x is defined at the middle of every bond 

along the original zigzag chain, x a = x + a — 3/2 define the sites of the subchains, MS is the ground 

state magnetization per spin, hp has a meaning of the Fermi momentum and is connected to M 

(e.g., kp = 7r(l + Af)/2 for S = 1/2), A; are nonuniversal amplitudes. All parameters (K, v, M, Ai) 

should be understood as functions of the field H which generally have to be extracted by comparing 

the numerical results with those following from the TLL description [ [17l [HI [HI [20l [2ll [23] . 

Introducing symmetric and antisymmetric combinations of the bosonic fields ip± — (ipi ± 

</32)/v / 2, 0± = (61 ±#2) / one obtains two coupled TLLs, and the nature of the leading interaction 

term strongly depends on the value of the TLL parameter K. The longitudinal (S Z S Z ) part of 

the zigzag exchange leads to a splitting of the TLL parameter values for the symmetric and 

1/2 

antisymmetric sectors, which for a ^> 1 are found as K± 1=3 A"[l ± 2K / (irva))\ , v± 1=3 v[l ± 

1/2 

2K/(TTva))] . At the same time, the S Z S Z part of the interaction gives rise to the interaction 
term proportional to cos ^/8n(p- — 7rM] whose scaling dimension is 2 AT_ and which favors n-type 
spin-nematic (or XY2 in the nomenclature of Ref. [[24]) correlations [[25]. The transversal part 



2 



Spin current in frustrated chains 



of the zigzag exchange yields the so-called "twist term" sin (^^/2n9-) (d x 8 + ) that has a nonzero 
conformal spin and the formal scaling dimension 1 + 1/(2KJ). The twist term favors states with 
spontaneous spin current [ [5]. If K < 1, as is the case for S = 1/2, the two interaction terms 
compete, which leads to the Ising-type transition between the ncmatic and spin-current state [ I12j . 
For S = 1 it is known that K > 1 [[5ni[Hl[23], thus the nematic-favoring term is irrelevant and the 
spin-current favoring term dominates. For higher spins, to our knowledge, there is no numerical 
data on the behavior of the TLL parameter K. However, since the analysis based on a mapping to 
the nonlinear a- model [|5D] can be extended to higher S, one may expect that for all spin-5 chains 
with S > 1, the TLL parameter K is greater than 1, at least in the region of small (but finite) 
magnetization M. 

Inside the phase with a spontaneous spin current, the twist term is the most relevant perturba- 
tion, and the mean-field treatment along the lines of Ref. [[9] leads to the conclusion [[12] that both 
(d x 0+) and (sin (\/27r#-)) should be nonzero, and the antisymmetric sector should be massive. 
Expressing the spin current k(ji) as defined in |T]) through the bosonic fields, one obtains [[5] 

k z {x) = sin (V2tt6-) Ia\ - {■KA l /2) 2 {d x 6 + ) 2 + (A 2 2 /2) cos (V8n(p + + 2k F x)}, 

n+(x) = 2A\M sin [^/ttJ29J) cxp {i[-kx/2 + y^/2 9 + ] }. (5) 

where we have omitted the contributions from massive fields and operators with higher scaling 
dimensions. At large distances x 3> £, where £ is the largest correlation length determined by the 
smallest gap in the neglected massive fields, the asymptotics of the leading spin-current correlations 
takes the form 

(k z (x)k z (0)) -> k§(1 + CiaT 4 + C 2 cos(2k F x)x- 4K+ ), 

(k + (x)k-(0)) -» AlM 2 x- 1/{iK+) exp{iQx}, (6) 

where the incommensurate wave vector Q is given by Q = tt/2 + ^Jn /2{d x 9 + ) and Cj are some 
numbers. 

The asymptotics ([6]), derived by means of bosonization, is expected to be valid in the limit 
a 3> 1 and for H sufficiently far from any of the critical fields H c , H s (close to H c or H s the 
bosonization approach becomes inapplicable since the effective bandwidth goes to zero and the 
interactions cease to be small compared to the bandwidth). In a close vicinity of the saturation 
field H s , large-iS analysis [[12] maps the system to a two-component dilute Bose gas with repulsive 
interaction, and the external magnetic field drives the system into a condensed state, playing the 
role of the chemical potential. The interspecies repulsion turns out to be strong enough to satisfy 
the phase separation condition, so only one of the species condenses and the other condensate is 
depleted. At the end, one deals with the one-component pseudo-condensate whose physics is again 
described by a (one-component) TLL, and the asymptotic form of the spin current correlators for 
H close to H s is given by [112): 

(k z (x)k z (0)) -^k 2 - C/x 2 , (k + (x)k-(0)) -» C'x- 1/{2K,) e lQ ' x . (7) 

where k\ cx (H s - H), C cx (H s - H) l / 2 - l /^ K '\ K' > 1 is another TLL parameter (different 
from K considered above), and the wave vector Q' is just the pitch of the classical helix, Q' = 
±(tt — arccos(l/4a)). 

3. Results of numerical study 

We have studied numerically the frustrated chain model ([2|) for systems with the spin S — 1/2, 
1, 3/2, and 2 and the frustration parameter a — J%j J\ was fixed at a = 1 (so the chains may be 
considered as stripes of the triangular lattice). We have used the density matrix renormalization 
group (DMRG) method [ 26, 27] in its matrix product state formulation [ 28J, making full use of the 
non-Abelian SU(2) symmetry [ [29J [28] . Although the external field breaks SU(2) symmetry, the 
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M=S t JLS 

tot 

Figure 1. Square of the ground state nearest-neighbor spin current «o (longitudinal vector chi- 
rality) as a function of the magnetization per spin M — Stot/LS, for isotropic spin-5 1 chains 
with the spin S = 1/2, 1, 3/2 and 2. The error bars shown correspond to the uncertainties of 
the fit. 



fact that the Zeeman energy term commutes with the rest of the Hamiltonian makes it possible 
to include the effect of the magnetic field by calculating the ground state in a sector with a 
finite total spin S to t- The use of SU(2) symmetry allows one to reduce substantially the number 
of states m which is necessary to describe the system: essentially, one treats the multiplet of 
states of the same total spin as a single representative state. A slight disadvantage is that the 
non-Abelian method allows to compute only reduced matrix elements (in the sense of the Wigner- 
Eckart theorem). Since the spin current is a vector, the available correlator that is easily available 
is the rotationally invariant scalar product ■ K,(n')), which is a mixture of the longitudinal 

and transversal spin current correlations. This complicates somewhat the analysis of the data: 
from the theoretical analysis in Sect. [Hit follows that for typical values of K ~ 1 the longitudinal 
correlations (k z (x)k z (0)} decay to their asymptotic value much faster than the transversal ones. 
Thus by analyzing the scalar correlator (/«(n) • k(ti')) one can only estimate the critical exponent 
r\ = T)x characterizing the decay of the transversal spin current correlations. 

We have studied zigzag chains of several lengths L ranging from 64 to 256. ft turned out that 
a relatively large number m of representative states was necessary to reach good convergence; 
depending on S, L and Stot, we have used m ranging from 400 to 1500. The scalar spin current 
correlator (re(n) • k(ti')) has been calculated for a large number of ground states in sectors with 
different total spin Stot- 

The correlator has been averaged over the starting and final positions n, n', and contributions 
with n or n' being closer as as a fixed "cutoff" (taken here to be 20 sites) to the chain ends were 
discarded. The DMRG data for the correlator has been fitted to the power-law form 

(k(x) ■ k(0)) = /eg + Ax ~ V cos[q(x + 5)] (8) 

suggested by (JB]), ©; the introduction of a finite phase shift 5 helps to suit the open boundary 
conditions. From those fits we have extracted the behavior of the equilibrium spin current kq and 
the exponent r] as functions of the chain magnetization M = Stot/L, shown respectively in Fig. [1] 
and Fig. [2] The wave vector q extracted from the fit depends only weakly on the magnetization, 
for all S studied, so we do not present the corresponding dependences. 

The quality of fits is generally at its best in the intermediate M region, and is deteriorating 
at M — * or M — > 1 for the reasons discussed in Ref. [ [B] . For small M it was only possible to 
extract the asymptotic value k,q of the spin current correlator, but not the critical exponent 77. 
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Figure 2. Critical exponent rj, which determines the power-law decay of the transversal spin 
current correlations, as a function of the magnetization per spin MM = Stat/LS, for isotropic 
spin-S chains with the spin S = 1/2, 1, 3/2 and 2. The error bars shown correspond to the 
uncertainties of the fit. 

The extracted value of r\ is weakly depending on M and for S > 1 qualitatively agrees with the 
theoretical estimates predicting that rj should vary from approximately 1/4 to 1/2 as M varies 
from to 1. One should realize that the error bars shown in Figs. ITffi are only of indicative nature, 
since they show just the uncertainties of the fit to the fixed form ((5} and do not take into account 
the variations of the fit parameters which would result from adding to it subleading contributions. 

4. Summary and discussion 

We have studied spin-S* isotropic antiferromagnetic frustrated chains in strong magnetic fields, 
described by the model |J2J) with the frustration parameter a = 1 and the spin S ranging from 1/2 
to 2, by means of the matrix product density matrix renormalization group technique. Existence 
of a phase with field-induced spontaneous spin current (vector chiral, or p-type nematic order) is 
established for all S studied. For S > 1, the behavior of the order parameter and its correlations 
as functions of the magnetization qualitatively agrees with the theory proposed in Ref. [ [12] . For 
S = 1/2, the existence of the vector chiral phase and the transition between the chiral and XY2 
phases are captured well by the large-cv bosonization approach of Ref. [ [12] . 

However, the behavior of S = 1/2 chain close to the saturation M = 1 does not fit into the 
common large- S scheme of [ 112]: there is a finite non-chiral region in the immediate vicinity of 
M = 1 which, as shown in [ 13J, belongs to the two-component Tomonaga-Luttinger liquid (TLL) 
phase. This peculiarity might be attributed to the fact that the large-iS approach in [[12] is based on 
identifying the two emerging particle species as bosons, while the underlying particles for S — 1/2 
are the Jordan- Wigner fermions. From the Bethe-ansatz results for the Hubbard model [ 30J it is 
known that pure density-density interaction does not destroy the two-component TLL in the low 
particle density limit. Thus, for S =1/2 the nature of the transition from the two-component TLL 
into the chiral phase close to the saturation field still remains to be explained. 

Acknowledgements 

We are grateful to T. Vekua for fruitful discussions. AK was supported by the Heisenberg 
Program (grant KO 2335/1-2) of Deutsche Forschungsgemcinschaft. 



5 



A.K. Kolezhuk, I. P. McCulloch 



References 

1. J. Villain, Ann. Isr. Phys. Soc. 2, 565 (1978). 

2. A.F.Andreev and LA. Grishchuk, Sov. Phys. JETP 60, 267 (1984). 

3. M. Kaburagi, H. Kawamura, and T. Hikihara, J. Phys. Soc. Jpn. 68, 3185 (1999). 

4. T. Hikihara, M. Kaburagi, and H. Kawamura, Phys. Rev. B 63, 174430 (2001). 

5. N. Shannon, T. Momoi, P. Sindzingre, Phys. Rev. Lett. 96, 027213 (2006). 

6. I. P. McCulloch, R. Kube, M. Kurz, A. Kleine, U. Schollwock, A. K. Kolezhuk, Phys. Rev. B 77, 
094404 (2008). 

7. T. Hikihara, L. Kecke, T. Momoi, A. Furusaki, Phys. Rev. B 78, 144404 (2008); J. Sudan, A. Luscher, 
A. Laeuchli, preprint larXiv:0807. 1923 1 

8. P. Chandra and P. Coleman, Phys. Rev. Lett. 66, 100 (1991). 

9. A. A. Nersesyan, A. O. Gogolin, and F. H. L. Efiler, Phys. Rev. Lett. 81, 910 (1998). 

10. A. K. Kolezhuk, Phys. Rev. B 62, R6057 (2000). 

11. P. Lecheminant, T. Jolicoeur, and P. Azaria, Phys. Rev. B 63, 174426 (2001). 

12. A. Kolezhuk and T. Vekua, Phys. Rev. B 72, 094424 (2005). 

13. K. Okunishi, J. Phys. Soc. Jpn. 77, 114004 (2008). 

14. M. Enderle, C. Mukherjee, B. Fak, R. K. Kremer, J.-M. Broto, H. Rosner, S.-L. Drechsler, J. Richter, 
J. Malek, A. Prokofiev, W. Assmus, S. Pujol, J.-L. Raggazzoni, H. Rakoto, M. Rheinstadter, and H. 
M. Ronnow, Europhys. Lett. 70, 237 (2005). 

15. S.-L. Drechsler, O. Volkova, A. N. Vasiliev, N. Tristan, J. Richter, M. Schmitt, H. Rosner, J. Malek, 
R. Klingeler, A. A. Zvyagin, and B. Biichner, Phys. Rev. Lett. 98, 077202 (2007). 

16. V. O. Garlea, A. Zheludev, K. Habicht, M. Meissner, B. Grenier, L.-P. Regnault, and E. Ressouche, 
Phys. Rev. B 79, 060404(R) (2009). 

17. N. M. Bogoliubov, A. G. Izergin, and V. E. Korepin, Nucl. Phys. B 275, 687 (1986). 

18. I. Affleck and M. Oshikawa, Phys. Rev. B 60, 1038 (1999). 

19. F. H. L. Essler, A. Furusaki, and T. Hikihara, Phys. Rev. B 68, 064410 (2003); T. Hikihara and A. 
Furusaki, ibid. 69, 064427 (2004). 

20. R. M. Konik and P. Fendley Phys. Rev. B 66, 144416 (2002). 

21. L. Campos Venuti, E. Ercolessi, G. Morandi, P. Pieri, and M. Roncaglia, Int. J. Mod. Phys. B 16, 
1363 (2002). 

22. Masahiro Sato, J. Stat. Mech. P09001 (2006). 

23. G. Fath, Phys. Rev. B 68, 134445 (2003). 

24. H. J. Schulz: Phys. Rev. B 34, 6372 (1986). 

25. T. Vekua, A. Honecker, H.-J. Mikeska, and F. Heidrich-Meisner, Phys. Rev. B 76, 174420 (2007). 

26. S. R. White, Phys. Rev. Lett. 69, 2863 (1992). 

27. U. Schollwock, Rev. Mod. Phys. 77, 259 (2005). 

28. I. P. McCulloch, J. Stat. Mech.: Theor. Exp., P10014 (2007). 

29. I. P. McCulloch and M. Gulacsi, Europhys. Lett. 57, 852 (2002). 

30. F. H. L. Essler, H. Frahm, F. Gohmann, A. Klumper, and V. E. Korepin, "The One-Dimensional 
Hubbard Model" (Cambridge University Press, 2005). 



6 



